Overview of analysis

This analysis examines if overall community composition differs based on self-reported ethnicity and what specific taxa are differentially abundant. MetaPhlAn (v3.0.7) was used to profile the taxonomic composition of metagenomes and calculate unweighted and weighted UniFrac distances.

Setting up the environment

Load the necessary libraries.

library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.0 ──
## ✓ ggplot2 3.3.2     ✓ purrr   0.3.4
## ✓ tibble  3.0.4     ✓ dplyr   1.0.2
## ✓ tidyr   1.1.2     ✓ stringr 1.4.0
## ✓ readr   1.4.0     ✓ forcats 0.5.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(vegan)
## Loading required package: permute
## Loading required package: lattice
## This is vegan 2.5-7
library(ggplot2)
library(cowplot)
library(car)
## Loading required package: carData
## 
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
## The following object is masked from 'package:purrr':
## 
##     some
library(multcomp)
## Loading required package: mvtnorm
## Loading required package: survival
## Loading required package: TH.data
## Loading required package: MASS
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
## 
## Attaching package: 'TH.data'
## The following object is masked from 'package:MASS':
## 
##     geyser
library(ggpubr)
## 
## Attaching package: 'ggpubr'
## The following object is masked from 'package:cowplot':
## 
##     get_legend

Importing data

Import metadata

metadata = read_tsv("year2_sample_mapping_metaphlan.tsv")
## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   .default = col_character(),
##   seq_id_simple = col_double(),
##   fecal_use = col_double(),
##   oral_use = col_double(),
##   mock = col_double(),
##   blank = col_double(),
##   day_number = col_double(),
##   time = col_time(format = ""),
##   extraction_kit_use = col_double(),
##   extraction = col_double(),
##   library_plate = col_double(),
##   id_num = col_double()
## )
## ℹ Use `spec()` for the full column specifications.

Import species-level relative abundance table, transpose, and recalculate relative abundance to sum to 1 instead of 100.

relab = read_tsv("year2_merged_abundance_table_species.txt")
## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   .default = col_double(),
##   body_site = col_character()
## )
## ℹ Use `spec()` for the full column specifications.
relab1 = relab %>% mutate_if(is.numeric, ~./100)
relab1_t = relab1 %>% gather(file, count, 2:164) %>% 
  spread(body_site, count)

Import UniFrac matrices

wuni_fecal = as.dist(
  read.table("year2_merged_abundance_wunifrac_fecal.tsv"))
uuni_fecal = as.dist(
  read.table("year2_merged_abundance_uunifrac_fecal.tsv"))
wuni_oral = as.dist(
  read.table("year2_merged_abundance_wunifrac_oral.tsv"))
uuni_oral = as.dist(
  read.table("year2_merged_abundance_uunifrac_oral.tsv"))

wuni_fecalb = as.dist(
  read.table("year2_merged_abundance_wunifrac_fecalb.tsv"))
uuni_fecalb = as.dist(
  read.table("year2_merged_abundance_uunifrac_fecalb.tsv"))
wuni_oralb = as.dist(
  read.table("year2_merged_abundance_wunifrac_oralb.tsv"))
uuni_oralb = as.dist(
  read.table("year2_merged_abundance_uunifrac_oralb.tsv"))

wuni_fecala = as.dist(
  read.table("year2_merged_abundance_wunifrac_fecala.tsv"))
uuni_fecala = as.dist(
  read.table("year2_merged_abundance_uunifrac_fecala.tsv"))
wuni_orala = as.dist(
  read.table("year2_merged_abundance_wunifrac_orala.tsv"))
uuni_orala = as.dist(
  read.table("year2_merged_abundance_uunifrac_orala.tsv"))

wuni_fecald = as.dist(
  read.table("year2_merged_abundance_wunifrac_fecald.tsv"))
uuni_fecald = as.dist(
  read.table("year2_merged_abundance_uunifrac_fecald.tsv"))

Filter samples

Filter out samples that have no identified species or are controls (blanks and mock communities). Separate oral samples from fecal samples.

relab_fecal = relab1_t %>% left_join(metadata, by = "file") %>% 
  filter(source == "fecal") 
relab_fecal_filtered = relab_fecal %>% 
  mutate(total = rowSums(relab_fecal[,2:496])) %>% 
  filter(total > 0)
relab_oral = relab1_t %>% left_join(metadata, by = "file") %>% 
  filter(source == "oral")

Calculating beta diversity metrics

UniFrac distances are calculated by a MetaPhlAn helper function. Vegan (2.5-7) is used to calculate the Bray-Curtis dissimilarity and Jaccard similarity metrics.

brayf = vegdist(relab_fecal_filtered[,2:496], method = "bray")
jaccf = vegdist(relab_fecal_filtered[,2:496], method = "jaccard")
brayo = vegdist(relab_oral[,2:496], method = "bray")
jacco = vegdist(relab_oral[,2:496], method = "jaccard")

PERMANOVAs - All Samples

Differences in the taxonomic structure of the community associated with self-reported ethnicity is assessed using permutational multivariate analysis of variance.

Bray-Curtis

Fecal

adonis2(brayf ~ ethnicity, data = relab_fecal_filtered, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayf ~ ethnicity, data = relab_fecal_filtered, permutations = 4999)
##            Df SumOfSqs      R2      F Pr(>F)    
## ethnicity   1   1.4525 0.06536 7.6225  2e-04 ***
## Residual  109  20.7702 0.93464                  
## Total     110  22.2227 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(brayo ~ ethnicity, data = relab_oral, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayo ~ ethnicity, data = relab_oral, permutations = 4999)
##           Df SumOfSqs     R2      F Pr(>F)   
## ethnicity  1   0.3486 0.0908 3.1957 0.0042 **
## Residual  32   3.4905 0.9092                 
## Total     33   3.8391 1.0000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Bray-Curtis dissimilarity with all samples.

## Run 0 stress 0.1730196 
## Run 1 stress 0.2313585 
## Run 2 stress 0.1900919 
## Run 3 stress 0.2471194 
## Run 4 stress 0.1900924 
## Run 5 stress 0.2006392 
## Run 6 stress 0.206551 
## Run 7 stress 0.1972424 
## Run 8 stress 0.2333067 
## Run 9 stress 0.2111622 
## Run 10 stress 0.1857487 
## Run 11 stress 0.2300488 
## Run 12 stress 0.173776 
## Run 13 stress 0.2264283 
## Run 14 stress 0.173776 
## Run 15 stress 0.192742 
## Run 16 stress 0.1744768 
## Run 17 stress 0.2105569 
## Run 18 stress 0.2032629 
## Run 19 stress 0.2129299 
## Run 20 stress 0.2036804 
## Run 21 stress 0.2029093 
## Run 22 stress 0.2086552 
## Run 23 stress 0.1867749 
## Run 24 stress 0.2201923 
## Run 25 stress 0.1777439 
## Run 26 stress 0.1730196 
## ... New best solution
## ... Procrustes: rmse 1.773753e-05  max resid 8.728979e-05 
## ... Similar to previous best
## *** Solution reached

## Run 0 stress 0.1592337 
## Run 1 stress 0.158507 
## ... New best solution
## ... Procrustes: rmse 0.009437815  max resid 0.04843542 
## Run 2 stress 0.158507 
## ... Procrustes: rmse 7.10493e-06  max resid 2.3685e-05 
## ... Similar to previous best
## Run 3 stress 0.1869209 
## Run 4 stress 0.3863172 
## Run 5 stress 0.1591203 
## Run 6 stress 0.1817134 
## Run 7 stress 0.1660616 
## Run 8 stress 0.1664842 
## Run 9 stress 0.1797591 
## Run 10 stress 0.1632234 
## Run 11 stress 0.1632234 
## Run 12 stress 0.1665592 
## Run 13 stress 0.1629265 
## Run 14 stress 0.1609917 
## Run 15 stress 0.1798961 
## Run 16 stress 0.1665592 
## Run 17 stress 0.2011488 
## Run 18 stress 0.1665592 
## Run 19 stress 0.1632234 
## Run 20 stress 0.179896 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_bray.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(brayf + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(brayf + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(brayo + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Jaccard

Fecal

adonis2(jaccf ~ ethnicity, data = relab_fecal_filtered, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccf ~ ethnicity, data = relab_fecal_filtered, permutations = 4999)
##            Df SumOfSqs      R2      F Pr(>F)    
## ethnicity   1    1.904 0.05969 6.9197  2e-04 ***
## Residual  109   29.992 0.94031                  
## Total     110   31.896 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(jacco ~ ethnicity, data = relab_oral, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jacco ~ ethnicity, data = relab_oral, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)   
## ethnicity  1   0.4911 0.07269 2.5084 0.0038 **
## Residual  32   6.2653 0.92731                 
## Total     33   6.7564 1.00000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Jaccard similarity with all samples.

## Run 0 stress 0.173776 
## Run 1 stress 0.1730198 
## ... New best solution
## ... Procrustes: rmse 0.03083177  max resid 0.1163365 
## Run 2 stress 0.2027132 
## Run 3 stress 0.1977463 
## Run 4 stress 0.2056424 
## Run 5 stress 0.2249779 
## Run 6 stress 0.2209721 
## Run 7 stress 0.2121922 
## Run 8 stress 0.2191132 
## Run 9 stress 0.2062637 
## Run 10 stress 0.1878283 
## Run 11 stress 0.2035525 
## Run 12 stress 0.2064486 
## Run 13 stress 0.2054211 
## Run 14 stress 0.2063329 
## Run 15 stress 0.1777443 
## Run 16 stress 0.2231832 
## Run 17 stress 0.2035501 
## Run 18 stress 0.2056861 
## Run 19 stress 0.1884853 
## Run 20 stress 0.1919388 
## Run 21 stress 0.1999187 
## Run 22 stress 0.2057736 
## Run 23 stress 0.173776 
## Run 24 stress 0.2056837 
## Run 25 stress 0.1973273 
## Run 26 stress 0.2022282 
## Run 27 stress 0.2091878 
## Run 28 stress 0.2174981 
## Run 29 stress 0.2326594 
## Run 30 stress 0.2082968 
## Run 31 stress 0.1974395 
## Run 32 stress 0.1950047 
## Run 33 stress 0.1827343 
## Run 34 stress 0.1930404 
## Run 35 stress 0.2053464 
## Run 36 stress 0.1880971 
## Run 37 stress 0.1857488 
## Run 38 stress 0.2166193 
## Run 39 stress 0.2192174 
## Run 40 stress 0.2357384 
## Run 41 stress 0.2166677 
## Run 42 stress 0.2171435 
## Run 43 stress 0.2279408 
## Run 44 stress 0.1730199 
## ... Procrustes: rmse 7.797651e-05  max resid 0.0003221037 
## ... Similar to previous best
## *** Solution reached

## Run 0 stress 0.1592336 
## Run 1 stress 0.1624179 
## Run 2 stress 0.160991 
## Run 3 stress 0.2056952 
## Run 4 stress 0.2032465 
## Run 5 stress 0.1835446 
## Run 6 stress 0.1609925 
## Run 7 stress 0.2177963 
## Run 8 stress 0.1930489 
## Run 9 stress 0.1629265 
## Run 10 stress 0.1665592 
## Run 11 stress 0.2176376 
## Run 12 stress 0.1622289 
## Run 13 stress 0.1968443 
## Run 14 stress 0.2073231 
## Run 15 stress 0.1823253 
## Run 16 stress 0.1676377 
## Run 17 stress 0.1629265 
## Run 18 stress 0.162229 
## Run 19 stress 0.175451 
## Run 20 stress 0.1754509 
## Run 21 stress 0.1676378 
## Run 22 stress 0.1591202 
## ... New best solution
## ... Procrustes: rmse 0.02981526  max resid 0.1456574 
## Run 23 stress 0.2009957 
## Run 24 stress 0.1609921 
## Run 25 stress 0.160991 
## Run 26 stress 0.1859658 
## Run 27 stress 0.1766694 
## Run 28 stress 0.1798963 
## Run 29 stress 0.1862181 
## Run 30 stress 0.1629265 
## Run 31 stress 0.1665592 
## Run 32 stress 0.1632234 
## Run 33 stress 0.2038537 
## Run 34 stress 0.162229 
## Run 35 stress 0.1632234 
## Run 36 stress 0.1835444 
## Run 37 stress 0.1609911 
## Run 38 stress 0.1629265 
## Run 39 stress 0.1622289 
## Run 40 stress 0.1609928 
## Run 41 stress 0.1699121 
## Run 42 stress 0.158507 
## ... New best solution
## ... Procrustes: rmse 0.03767444  max resid 0.1911456 
## Run 43 stress 0.2017916 
## Run 44 stress 0.2106544 
## Run 45 stress 0.158507 
## ... New best solution
## ... Procrustes: rmse 1.27024e-05  max resid 4.391212e-05 
## ... Similar to previous best
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_jacc.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(jaccf + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(jaccf + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(jacco + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Weighted UniFrac

Fecal

adonis2(wuni_fecal ~ ethnicity, data = relab_fecal_filtered, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_fecal ~ ethnicity, data = relab_fecal_filtered, permutations = 4999)
##            Df SumOfSqs      R2      F Pr(>F)   
## ethnicity   1   0.2429 0.04942 5.6666 0.0012 **
## Residual  109   4.6721 0.95058                 
## Total     110   4.9150 1.00000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(wuni_oral ~ ethnicity, data = relab_oral, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_oral ~ ethnicity, data = relab_oral, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.07415 0.04607 1.5456 0.1934
## Residual  32  1.53533 0.95393              
## Total     33  1.60948 1.00000

NMDS visualization of weighted Unifrac distances with all samples.

## Run 0 stress 0.1128885 
## Run 1 stress 0.1210185 
## Run 2 stress 0.1144106 
## Run 3 stress 0.1144077 
## Run 4 stress 0.130739 
## Run 5 stress 0.1129022 
## ... Procrustes: rmse 0.0009940589  max resid 0.006519917 
## ... Similar to previous best
## Run 6 stress 0.1168925 
## Run 7 stress 0.1128926 
## ... Procrustes: rmse 0.006453419  max resid 0.0635843 
## Run 8 stress 0.1128787 
## ... New best solution
## ... Procrustes: rmse 0.001119216  max resid 0.009342401 
## ... Similar to previous best
## Run 9 stress 0.1163883 
## Run 10 stress 0.1380249 
## Run 11 stress 0.1486133 
## Run 12 stress 0.1425132 
## Run 13 stress 0.1129707 
## ... Procrustes: rmse 0.006421709  max resid 0.06307331 
## Run 14 stress 0.1352977 
## Run 15 stress 0.113063 
## ... Procrustes: rmse 0.009634561  max resid 0.07258814 
## Run 16 stress 0.1132554 
## ... Procrustes: rmse 0.00704866  max resid 0.07112905 
## Run 17 stress 0.1144078 
## Run 18 stress 0.1129007 
## ... Procrustes: rmse 0.001398501  max resid 0.009298345 
## ... Similar to previous best
## Run 19 stress 0.1204263 
## Run 20 stress 0.1128765 
## ... New best solution
## ... Procrustes: rmse 0.0007308774  max resid 0.005945148 
## ... Similar to previous best
## *** Solution reached

## Run 0 stress 0.08960119 
## Run 1 stress 0.09882517 
## Run 2 stress 0.1031962 
## Run 3 stress 0.1126862 
## Run 4 stress 0.115512 
## Run 5 stress 0.1086049 
## Run 6 stress 0.1073695 
## Run 7 stress 0.1027136 
## Run 8 stress 0.1179109 
## Run 9 stress 0.1077257 
## Run 10 stress 0.0896885 
## ... Procrustes: rmse 0.01592665  max resid 0.07752264 
## Run 11 stress 0.08960119 
## ... Procrustes: rmse 2.834432e-05  max resid 0.0001151919 
## ... Similar to previous best
## Run 12 stress 0.0896885 
## ... Procrustes: rmse 0.01590753  max resid 0.07742581 
## Run 13 stress 0.1027136 
## Run 14 stress 0.1030229 
## Run 15 stress 0.09882528 
## Run 16 stress 0.1025333 
## Run 17 stress 0.1086048 
## Run 18 stress 0.08968854 
## ... Procrustes: rmse 0.01593447  max resid 0.0774934 
## Run 19 stress 0.1004283 
## Run 20 stress 0.1085053 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_wuni.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(wunif + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(wunif + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(wunio + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Unweighted UniFrac

Fecal

adonis2(uuni_fecal ~ ethnicity, data = relab_fecal_filtered, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_fecal ~ ethnicity, data = relab_fecal_filtered, permutations = 4999)
##            Df SumOfSqs      R2      F Pr(>F)    
## ethnicity   1   0.7732 0.07238 8.5048  2e-04 ***
## Residual  109   9.9100 0.92762                  
## Total     110  10.6832 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(uuni_oral ~ ethnicity, data = relab_oral, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_oral ~ ethnicity, data = relab_oral, permutations = 4999)
##           Df SumOfSqs     R2      F Pr(>F)
## ethnicity  1  0.05327 0.0402 1.3401 0.1886
## Residual  32  1.27193 0.9598              
## Total     33  1.32520 1.0000

NMDS visualization of Bray-Curtis dissimilarity with all samples.

## Run 0 stress 7.26794e-05 
## Run 1 stress 9.69908e-05 
## ... Procrustes: rmse 3.549433e-05  max resid 0.0001432463 
## ... Similar to previous best
## Run 2 stress 5.901279e-05 
## ... New best solution
## ... Procrustes: rmse 3.197825e-05  max resid 9.467394e-05 
## ... Similar to previous best
## Run 3 stress 9.395146e-05 
## ... Procrustes: rmse 4.761152e-05  max resid 0.0001204205 
## ... Similar to previous best
## Run 4 stress 9.54084e-05 
## ... Procrustes: rmse 2.881724e-05  max resid 8.458948e-05 
## ... Similar to previous best
## Run 5 stress 7.52498e-05 
## ... Procrustes: rmse 2.551549e-05  max resid 7.243461e-05 
## ... Similar to previous best
## Run 6 stress 7.467214e-05 
## ... Procrustes: rmse 2.395832e-05  max resid 5.006309e-05 
## ... Similar to previous best
## Run 7 stress 9.933767e-05 
## ... Procrustes: rmse 3.384435e-05  max resid 0.0001055956 
## ... Similar to previous best
## Run 8 stress 6.422304e-05 
## ... Procrustes: rmse 2.370323e-05  max resid 6.051647e-05 
## ... Similar to previous best
## Run 9 stress 6.900282e-05 
## ... Procrustes: rmse 2.039856e-05  max resid 4.60642e-05 
## ... Similar to previous best
## Run 10 stress 9.428947e-05 
## ... Procrustes: rmse 2.37946e-05  max resid 6.96226e-05 
## ... Similar to previous best
## Run 11 stress 9.90301e-05 
## ... Procrustes: rmse 2.512668e-05  max resid 7.256941e-05 
## ... Similar to previous best
## Run 12 stress 8.85085e-05 
## ... Procrustes: rmse 2.767042e-05  max resid 7.239274e-05 
## ... Similar to previous best
## Run 13 stress 7.884904e-05 
## ... Procrustes: rmse 1.963874e-05  max resid 6.988197e-05 
## ... Similar to previous best
## Run 14 stress 9.487501e-05 
## ... Procrustes: rmse 2.634122e-05  max resid 7.477966e-05 
## ... Similar to previous best
## Run 15 stress 9.233808e-05 
## ... Procrustes: rmse 3.330246e-05  max resid 9.539123e-05 
## ... Similar to previous best
## Run 16 stress 9.531369e-05 
## ... Procrustes: rmse 3.363245e-05  max resid 0.000101479 
## ... Similar to previous best
## Run 17 stress 9.688617e-05 
## ... Procrustes: rmse 2.839829e-05  max resid 6.846969e-05 
## ... Similar to previous best
## Run 18 stress 9.723656e-05 
## ... Procrustes: rmse 4.005544e-05  max resid 8.558437e-05 
## ... Similar to previous best
## Run 19 stress 8.090206e-05 
## ... Procrustes: rmse 2.143516e-05  max resid 6.771863e-05 
## ... Similar to previous best
## Run 20 stress 8.926949e-05 
## ... Procrustes: rmse 3.082847e-05  max resid 9.812705e-05 
## ... Similar to previous best
## *** Solution reached
## Warning in metaMDS(uuni_fecal, k = 2, trymax = 500): stress is (nearly) zero:
## you may have insufficient data
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Run 0 stress 0.1791554 
## Run 1 stress 0.1822924 
## Run 2 stress 0.1954749 
## Run 3 stress 0.1806555 
## Run 4 stress 0.1839397 
## Run 5 stress 0.1811046 
## Run 6 stress 0.1817662 
## Run 7 stress 0.2109138 
## Run 8 stress 0.1829339 
## Run 9 stress 0.1791029 
## ... New best solution
## ... Procrustes: rmse 0.007979483  max resid 0.03643087 
## Run 10 stress 0.1839398 
## Run 11 stress 0.1792178 
## ... Procrustes: rmse 0.0194101  max resid 0.05380965 
## Run 12 stress 0.1822924 
## Run 13 stress 0.1825586 
## Run 14 stress 0.1933016 
## Run 15 stress 0.1791553 
## ... Procrustes: rmse 0.007764313  max resid 0.03549337 
## Run 16 stress 0.1809487 
## Run 17 stress 0.186249 
## Run 18 stress 0.1912585 
## Run 19 stress 0.1954286 
## Run 20 stress 0.1833687 
## Run 21 stress 0.1829339 
## Run 22 stress 0.1826122 
## Run 23 stress 0.1810658 
## Run 24 stress 0.1792178 
## ... Procrustes: rmse 0.0193816  max resid 0.05387272 
## Run 25 stress 0.1814989 
## Run 26 stress 0.1819031 
## Run 27 stress 0.1933562 
## Run 28 stress 0.182295 
## Run 29 stress 0.1791553 
## ... Procrustes: rmse 0.007764574  max resid 0.03552299 
## Run 30 stress 0.1787822 
## ... New best solution
## ... Procrustes: rmse 0.01530204  max resid 0.05363203 
## Run 31 stress 0.2137706 
## Run 32 stress 0.2208541 
## Run 33 stress 0.1787829 
## ... Procrustes: rmse 0.0007134825  max resid 0.002664038 
## ... Similar to previous best
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_uuni.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(uunif + theme(legend.box.margin = margin(0, 0, 0, 12)))
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
col1 = plot_grid(uunif + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
col2 = plot_grid(uunio + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

PERMANOVAs - Before samples

Filtering samples

Filter samples from before the diet intervention and recalculate Bray-Curtis and Jaccard metrics

relab_fecal_filtered_before = relab_fecal_filtered %>% 
  filter(Study_period == "Before")
relab_oral_before = relab_oral %>% 
  filter(Study_period == "Before")

brayfb = vegdist(relab_fecal_filtered_before[,2:496], method = "bray")
jaccfb = vegdist(relab_fecal_filtered_before[,2:496], method = "jaccard")
brayob = vegdist(relab_oral_before[,2:496], method = "bray")
jaccob = vegdist(relab_oral_before[,2:496], method = "jaccard")

Bray-Curtis

Fecal

adonis2(brayfb ~ ethnicity, data = relab_fecal_filtered_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayfb ~ ethnicity, data = relab_fecal_filtered_before, permutations = 4999)
##           Df SumOfSqs      R2     F Pr(>F)
## ethnicity  1   0.2744 0.06875 1.255  0.221
## Residual  17   3.7171 0.93125             
## Total     18   3.9915 1.00000

Oral

adonis2(brayob ~ ethnicity, data = relab_oral_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayob ~ ethnicity, data = relab_oral_before, permutations = 4999)
##           Df SumOfSqs      R2     F Pr(>F)
## ethnicity  1  0.18873 0.09118 1.505 0.1206
## Residual  15  1.88112 0.90882             
## Total     16  2.06985 1.00000

NMDS visualization of Bray-Curtis dissimilarity with before samples.

## Run 0 stress 0.1284112 
## Run 1 stress 0.1689742 
## Run 2 stress 0.1284112 
## ... Procrustes: rmse 1.267885e-06  max resid 2.283503e-06 
## ... Similar to previous best
## Run 3 stress 0.1871944 
## Run 4 stress 0.1697207 
## Run 5 stress 0.1906031 
## Run 6 stress 0.2872607 
## Run 7 stress 0.1525329 
## Run 8 stress 0.1284112 
## ... Procrustes: rmse 3.129066e-06  max resid 6.608683e-06 
## ... Similar to previous best
## Run 9 stress 0.1792602 
## Run 10 stress 0.1522154 
## Run 11 stress 0.1689742 
## Run 12 stress 0.1284112 
## ... Procrustes: rmse 1.700997e-05  max resid 3.624782e-05 
## ... Similar to previous best
## Run 13 stress 0.1525329 
## Run 14 stress 0.1792602 
## Run 15 stress 0.1444672 
## Run 16 stress 0.2020777 
## Run 17 stress 0.1444672 
## Run 18 stress 0.1284112 
## ... New best solution
## ... Procrustes: rmse 1.331546e-06  max resid 2.666263e-06 
## ... Similar to previous best
## Run 19 stress 0.1525331 
## Run 20 stress 0.1506507 
## *** Solution reached

## Run 0 stress 0.1253357 
## Run 1 stress 0.1300761 
## Run 2 stress 0.1384922 
## Run 3 stress 0.1364623 
## Run 4 stress 0.1264877 
## Run 5 stress 0.1467916 
## Run 6 stress 0.1294473 
## Run 7 stress 0.1300761 
## Run 8 stress 0.1453194 
## Run 9 stress 0.1453195 
## Run 10 stress 0.1264877 
## Run 11 stress 0.1461133 
## Run 12 stress 0.1256167 
## ... Procrustes: rmse 0.02122912  max resid 0.0615593 
## Run 13 stress 0.1300761 
## Run 14 stress 0.1295583 
## Run 15 stress 0.1253358 
## ... Procrustes: rmse 0.0002664434  max resid 0.000872116 
## ... Similar to previous best
## Run 16 stress 0.1455176 
## Run 17 stress 0.1264877 
## Run 18 stress 0.1384926 
## Run 19 stress 0.1295583 
## Run 20 stress 0.1364623 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_bray_before.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(brayfb + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(brayfb + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(brayob + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Jaccard

Fecal

adonis2(jaccfb ~ ethnicity, data = relab_fecal_filtered_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccfb ~ ethnicity, data = relab_fecal_filtered_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1   0.4122 0.07421 1.3628 0.1418
## Residual  17   5.1416 0.92579              
## Total     18   5.5537 1.00000

Oral

adonis2(jaccob ~ ethnicity, data = relab_oral_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccob ~ ethnicity, data = relab_oral_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1   0.2815 0.07968 1.2988 0.1374
## Residual  15   3.2508 0.92032              
## Total     16   3.5323 1.00000

NMDS visualization of Jaccard similarity with before samples.

## Run 0 stress 0.1284112 
## Run 1 stress 0.1284112 
## ... Procrustes: rmse 5.026012e-06  max resid 1.072812e-05 
## ... Similar to previous best
## Run 2 stress 0.1284112 
## ... Procrustes: rmse 2.255768e-06  max resid 4.757595e-06 
## ... Similar to previous best
## Run 3 stress 0.1643886 
## Run 4 stress 0.1689742 
## Run 5 stress 0.1444672 
## Run 6 stress 0.1284112 
## ... Procrustes: rmse 3.127735e-06  max resid 6.192814e-06 
## ... Similar to previous best
## Run 7 stress 0.1525329 
## Run 8 stress 0.1284112 
## ... Procrustes: rmse 8.96628e-07  max resid 1.623914e-06 
## ... Similar to previous best
## Run 9 stress 0.194464 
## Run 10 stress 0.152533 
## Run 11 stress 0.1697207 
## Run 12 stress 0.1284112 
## ... New best solution
## ... Procrustes: rmse 4.819807e-07  max resid 9.644547e-07 
## ... Similar to previous best
## Run 13 stress 0.1610714 
## Run 14 stress 0.1842235 
## Run 15 stress 0.2633195 
## Run 16 stress 0.1842236 
## Run 17 stress 0.1444672 
## Run 18 stress 0.152533 
## Run 19 stress 0.1689742 
## Run 20 stress 0.1945795 
## *** Solution reached

## Run 0 stress 0.1253358 
## Run 1 stress 0.1481009 
## Run 2 stress 0.1498922 
## Run 3 stress 0.1585965 
## Run 4 stress 0.1453194 
## Run 5 stress 0.1463235 
## Run 6 stress 0.1384927 
## Run 7 stress 0.1455176 
## Run 8 stress 0.1264877 
## Run 9 stress 0.1463234 
## Run 10 stress 0.1461133 
## Run 11 stress 0.1467916 
## Run 12 stress 0.1455177 
## Run 13 stress 0.1453195 
## Run 14 stress 0.1467915 
## Run 15 stress 0.1253357 
## ... New best solution
## ... Procrustes: rmse 0.0001619306  max resid 0.0005307169 
## ... Similar to previous best
## Run 16 stress 0.1453195 
## Run 17 stress 0.1253357 
## ... Procrustes: rmse 0.000132518  max resid 0.000433773 
## ... Similar to previous best
## Run 18 stress 0.1467916 
## Run 19 stress 0.1461134 
## Run 20 stress 0.1453194 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_jacc_before.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(jaccfb + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(jaccfb + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(jaccob + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Weighted UniFrac

Fecal

adonis2(wuni_fecalb ~ ethnicity, data = relab_fecal_filtered_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_fecalb ~ ethnicity, data = relab_fecal_filtered_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.02942 0.03414 0.6009 0.6462
## Residual  17  0.83220 0.96586              
## Total     18  0.86162 1.00000

Oral

adonis2(wuni_oralb ~ ethnicity, data = relab_oral_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_oralb ~ ethnicity, data = relab_oral_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.04241 0.04763 0.7501 0.4672
## Residual  15  0.84808 0.95237              
## Total     16  0.89049 1.00000

NMDS visualization of weighted Unifrac distances with before samples.

## Run 0 stress 0.07991057 
## Run 1 stress 0.08273946 
## Run 2 stress 0.08273945 
## Run 3 stress 0.07999749 
## ... Procrustes: rmse 0.009259724  max resid 0.02902095 
## Run 4 stress 0.08403418 
## Run 5 stress 0.08148165 
## Run 6 stress 0.08148163 
## Run 7 stress 0.08403424 
## Run 8 stress 0.08273944 
## Run 9 stress 0.07991038 
## ... New best solution
## ... Procrustes: rmse 0.000109443  max resid 0.0003392373 
## ... Similar to previous best
## Run 10 stress 0.08403433 
## Run 11 stress 0.07999758 
## ... Procrustes: rmse 0.009292403  max resid 0.02906161 
## Run 12 stress 0.08148165 
## Run 13 stress 0.2041272 
## Run 14 stress 0.1479628 
## Run 15 stress 0.1419221 
## Run 16 stress 0.2470987 
## Run 17 stress 0.2321205 
## Run 18 stress 0.08403426 
## Run 19 stress 0.08273946 
## Run 20 stress 0.08148165 
## *** Solution reached

## Run 0 stress 0.06894997 
## Run 1 stress 0.06966296 
## Run 2 stress 0.06894997 
## ... Procrustes: rmse 9.138811e-06  max resid 2.108204e-05 
## ... Similar to previous best
## Run 3 stress 0.07415722 
## Run 4 stress 0.07643535 
## Run 5 stress 0.07415714 
## Run 6 stress 0.07616554 
## Run 7 stress 0.07415714 
## Run 8 stress 0.07415716 
## Run 9 stress 0.07615363 
## Run 10 stress 0.07415714 
## Run 11 stress 0.06894999 
## ... Procrustes: rmse 2.244914e-05  max resid 5.319055e-05 
## ... Similar to previous best
## Run 12 stress 0.07415714 
## Run 13 stress 0.06894997 
## ... Procrustes: rmse 5.549988e-05  max resid 0.0001784903 
## ... Similar to previous best
## Run 14 stress 0.06894997 
## ... Procrustes: rmse 1.116434e-05  max resid 3.268703e-05 
## ... Similar to previous best
## Run 15 stress 0.07415714 
## Run 16 stress 0.06894997 
## ... Procrustes: rmse 6.667681e-06  max resid 1.536305e-05 
## ... Similar to previous best
## Run 17 stress 0.06966296 
## Run 18 stress 0.07047077 
## Run 19 stress 0.06895003 
## ... Procrustes: rmse 0.0001233991  max resid 0.0003953277 
## ... Similar to previous best
## Run 20 stress 0.06894997 
## ... New best solution
## ... Procrustes: rmse 2.332943e-05  max resid 7.447102e-05 
## ... Similar to previous best
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_wuni_before.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(wunifb + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(wunifb + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(wuniob + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Unweighted UniFrac

Fecal

adonis2(uuni_fecalb ~ ethnicity, data = relab_fecal_filtered_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_fecalb ~ ethnicity, data = relab_fecal_filtered_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1  0.18751 0.08793 1.6389 0.0192 *
## Residual  17  1.94502 0.91207                
## Total     18  2.13253 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(uuni_oralb ~ ethnicity, data = relab_oral_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_oralb ~ ethnicity, data = relab_oral_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.03204 0.04827 0.7608 0.7082
## Residual  15  0.63162 0.95173              
## Total     16  0.66365 1.00000

NMDS visualization of unweighted UniFrac distances with before samples.

## Run 0 stress 9.874988e-05 
## Run 1 stress 9.738051e-05 
## ... New best solution
## ... Procrustes: rmse 5.287388e-05  max resid 0.0001821857 
## ... Similar to previous best
## Run 2 stress 9.567545e-05 
## ... New best solution
## ... Procrustes: rmse 0.0001116539  max resid 0.0002738316 
## ... Similar to previous best
## Run 3 stress 9.531762e-05 
## ... New best solution
## ... Procrustes: rmse 0.0001555189  max resid 0.0003845121 
## ... Similar to previous best
## Run 4 stress 9.671068e-05 
## ... Procrustes: rmse 0.0001344295  max resid 0.0003474037 
## ... Similar to previous best
## Run 5 stress 5.434395e-05 
## ... New best solution
## ... Procrustes: rmse 0.000100357  max resid 0.0002018783 
## ... Similar to previous best
## Run 6 stress 9.900469e-05 
## ... Procrustes: rmse 0.0001047454  max resid 0.0001988247 
## ... Similar to previous best
## Run 7 stress 9.291372e-05 
## ... Procrustes: rmse 9.765583e-05  max resid 0.0002577777 
## ... Similar to previous best
## Run 8 stress 8.722854e-05 
## ... Procrustes: rmse 7.334147e-05  max resid 0.0001892078 
## ... Similar to previous best
## Run 9 stress 7.584566e-05 
## ... Procrustes: rmse 5.919929e-05  max resid 0.0001497131 
## ... Similar to previous best
## Run 10 stress 9.377897e-05 
## ... Procrustes: rmse 0.0001069664  max resid 0.0002615232 
## ... Similar to previous best
## Run 11 stress 9.360577e-05 
## ... Procrustes: rmse 0.0001182038  max resid 0.0003032233 
## ... Similar to previous best
## Run 12 stress 9.431176e-05 
## ... Procrustes: rmse 7.99195e-05  max resid 0.0001745545 
## ... Similar to previous best
## Run 13 stress 9.616363e-05 
## ... Procrustes: rmse 0.0001205418  max resid 0.0003031172 
## ... Similar to previous best
## Run 14 stress 6.990469e-05 
## ... Procrustes: rmse 5.261081e-05  max resid 0.0001056664 
## ... Similar to previous best
## Run 15 stress 9.423776e-05 
## ... Procrustes: rmse 9.508213e-05  max resid 0.0001781795 
## ... Similar to previous best
## Run 16 stress 9.876871e-05 
## ... Procrustes: rmse 8.801927e-05  max resid 0.0002047066 
## ... Similar to previous best
## Run 17 stress 9.25114e-05 
## ... Procrustes: rmse 7.214541e-05  max resid 0.0001535098 
## ... Similar to previous best
## Run 18 stress 9.256526e-05 
## ... Procrustes: rmse 0.0001124195  max resid 0.000271437 
## ... Similar to previous best
## Run 19 stress 8.026396e-05 
## ... Procrustes: rmse 7.486417e-05  max resid 0.0001550846 
## ... Similar to previous best
## Run 20 stress 8.552818e-05 
## ... Procrustes: rmse 9.18725e-05  max resid 0.0002441121 
## ... Similar to previous best
## *** Solution reached
## Warning in metaMDS(uuni_fecalb, k = 2, trymax = 500): stress is (nearly) zero:
## you may have insufficient data
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Run 0 stress 0.1694422 
## Run 1 stress 0.1715766 
## Run 2 stress 0.1703132 
## Run 3 stress 0.2433225 
## Run 4 stress 0.1694423 
## ... Procrustes: rmse 6.586859e-05  max resid 0.0001790191 
## ... Similar to previous best
## Run 5 stress 0.1715766 
## Run 6 stress 0.1694422 
## ... Procrustes: rmse 1.201997e-05  max resid 3.246919e-05 
## ... Similar to previous best
## Run 7 stress 0.2186982 
## Run 8 stress 0.1768019 
## Run 9 stress 0.1715766 
## Run 10 stress 0.1694423 
## ... Procrustes: rmse 4.487855e-05  max resid 0.0001187343 
## ... Similar to previous best
## Run 11 stress 0.2136743 
## Run 12 stress 0.1768019 
## Run 13 stress 0.1694423 
## ... Procrustes: rmse 3.643282e-05  max resid 0.0001001615 
## ... Similar to previous best
## Run 14 stress 0.1715766 
## Run 15 stress 0.1696274 
## ... Procrustes: rmse 0.01855263  max resid 0.05428006 
## Run 16 stress 0.1694423 
## ... Procrustes: rmse 0.0001232676  max resid 0.0003308309 
## ... Similar to previous best
## Run 17 stress 0.1715766 
## Run 18 stress 0.2199095 
## Run 19 stress 0.1715766 
## Run 20 stress 0.2198418 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_uuni_before.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(uunifb + theme(legend.box.margin = margin(0, 0, 0, 12)))
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
col1 = plot_grid(uunifb + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
col2 = plot_grid(uuniob + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

PERMANOVAs - During samples

Filter samples

Filter samples from during the diet intervention and recalculate Bray-Curtis and Jaccard metrics

relab_fecal_filtered_during = relab_fecal_filtered %>% 
  filter(Study_period == "During")

brayfd = vegdist(relab_fecal_filtered_during[,2:496], method = "bray")
jaccfd = vegdist(relab_fecal_filtered_during[,2:496], method = "jaccard")

Bray-Curtis

Fecal

adonis2(brayfd ~ ethnicity, data = relab_fecal_filtered_during, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayfd ~ ethnicity, data = relab_fecal_filtered_during, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)    
## ethnicity  1   0.7899 0.06703 4.2387  4e-04 ***
## Residual  59  10.9957 0.93297                  
## Total     60  11.7856 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Bray-Curtis dissimilarity with during samples.

## Run 0 stress 0.1664788 
## Run 1 stress 0.1752781 
## Run 2 stress 0.2053094 
## Run 3 stress 0.1998895 
## Run 4 stress 0.1752755 
## Run 5 stress 0.1969035 
## Run 6 stress 0.2004216 
## Run 7 stress 0.2313365 
## Run 8 stress 0.1981598 
## Run 9 stress 0.2247829 
## Run 10 stress 0.2112386 
## Run 11 stress 0.217677 
## Run 12 stress 0.224643 
## Run 13 stress 0.2069444 
## Run 14 stress 0.2063212 
## Run 15 stress 0.1980561 
## Run 16 stress 0.1958214 
## Run 17 stress 0.1664617 
## ... New best solution
## ... Procrustes: rmse 0.001075342  max resid 0.006206652 
## ... Similar to previous best
## Run 18 stress 0.1982036 
## Run 19 stress 0.2112386 
## Run 20 stress 0.1724189 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_bray_during.tif", res = 150, width = 11, 
     height = 7, units = "in")
brayfd
dev.off()
## quartz_off_screen 
##                 2

Jaccard

Fecal

adonis2(jaccfd ~ ethnicity, data = relab_fecal_filtered_during, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccfd ~ ethnicity, data = relab_fecal_filtered_during, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)    
## ethnicity  1   1.0741 0.06259 3.9392  2e-04 ***
## Residual  59  16.0879 0.93741                  
## Total     60  17.1620 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Jaccard similarity with during samples.

## Run 0 stress 0.1664788 
## Run 1 stress 0.2118375 
## Run 2 stress 0.1724188 
## Run 3 stress 0.175278 
## Run 4 stress 0.2004768 
## Run 5 stress 0.1664617 
## ... New best solution
## ... Procrustes: rmse 0.00107627  max resid 0.006162014 
## ... Similar to previous best
## Run 6 stress 0.1970234 
## Run 7 stress 0.1664788 
## ... Procrustes: rmse 0.0010852  max resid 0.006259243 
## ... Similar to previous best
## Run 8 stress 0.2210808 
## Run 9 stress 0.20746 
## Run 10 stress 0.1806977 
## Run 11 stress 0.1664617 
## ... New best solution
## ... Procrustes: rmse 3.550179e-05  max resid 0.0001502291 
## ... Similar to previous best
## Run 12 stress 0.2102412 
## Run 13 stress 0.1664617 
## ... Procrustes: rmse 1.514577e-05  max resid 4.866654e-05 
## ... Similar to previous best
## Run 14 stress 0.2057964 
## Run 15 stress 0.1998369 
## Run 16 stress 0.2055727 
## Run 17 stress 0.2008209 
## Run 18 stress 0.175278 
## Run 19 stress 0.2062881 
## Run 20 stress 0.1945864 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_jacc_during.tif", res = 150, width = 11, 
     height = 7, units = "in")
jaccfd
dev.off()
## quartz_off_screen 
##                 2

Weighted UniFrac

Fecal

adonis2(wuni_fecald ~ ethnicity, data = relab_fecal_filtered_during, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_fecald ~ ethnicity, data = relab_fecal_filtered_during, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)   
## ethnicity  1  0.15963 0.06033 3.7881  0.009 **
## Residual  59  2.48636 0.93967                 
## Total     60  2.64599 1.00000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of weighted Unifrac distances with during samples.

## Run 0 stress 0.1065766 
## Run 1 stress 0.1066787 
## ... Procrustes: rmse 0.01085342  max resid 0.08224327 
## Run 2 stress 0.1506047 
## Run 3 stress 0.1159546 
## Run 4 stress 0.1167693 
## Run 5 stress 0.1164082 
## Run 6 stress 0.1412151 
## Run 7 stress 0.116769 
## Run 8 stress 0.1168791 
## Run 9 stress 0.1168802 
## Run 10 stress 0.143699 
## Run 11 stress 0.1182275 
## Run 12 stress 0.1157516 
## Run 13 stress 0.1363755 
## Run 14 stress 0.1182267 
## Run 15 stress 0.1405206 
## Run 16 stress 0.1125866 
## Run 17 stress 0.1103593 
## Run 18 stress 0.1459305 
## Run 19 stress 0.1186453 
## Run 20 stress 0.1371618 
## Run 21 stress 0.1103092 
## Run 22 stress 0.1103094 
## Run 23 stress 0.1066172 
## ... Procrustes: rmse 0.004435014  max resid 0.02910857 
## Run 24 stress 0.1066786 
## ... Procrustes: rmse 0.01084442  max resid 0.08216484 
## Run 25 stress 0.1159532 
## Run 26 stress 0.1170446 
## Run 27 stress 0.1404304 
## Run 28 stress 0.1159534 
## Run 29 stress 0.1066787 
## ... Procrustes: rmse 0.0108583  max resid 0.08228771 
## Run 30 stress 0.1339132 
## Run 31 stress 0.1295888 
## Run 32 stress 0.1170402 
## Run 33 stress 0.1066172 
## ... Procrustes: rmse 0.004433051  max resid 0.02912203 
## Run 34 stress 0.1065764 
## ... New best solution
## ... Procrustes: rmse 8.767366e-05  max resid 0.0004006176 
## ... Similar to previous best
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_wuni_during.tif", res = 150, width = 11, 
     height = 7, units = "in")
wunif
dev.off()
## quartz_off_screen 
##                 2

Unweighted UniFrac

Fecal

adonis2(uuni_fecald ~ ethnicity, data = relab_fecal_filtered_during, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_fecald ~ ethnicity, data = relab_fecal_filtered_during, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)    
## ethnicity  1   0.4436 0.07779 4.9764  2e-04 ***
## Residual  59   5.2595 0.92221                  
## Total     60   5.7031 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of unweighted Unifrac distance with during samples.

## Run 0 stress 0.1793941 
## Run 1 stress 0.2142523 
## Run 2 stress 0.2151973 
## Run 3 stress 0.1848087 
## Run 4 stress 0.1793941 
## ... New best solution
## ... Procrustes: rmse 7.947352e-05  max resid 0.0004143951 
## ... Similar to previous best
## Run 5 stress 0.1847541 
## Run 6 stress 0.2180861 
## Run 7 stress 0.1815689 
## Run 8 stress 0.1794329 
## ... Procrustes: rmse 0.006333638  max resid 0.03047401 
## Run 9 stress 0.1862549 
## Run 10 stress 0.1800931 
## Run 11 stress 0.1869278 
## Run 12 stress 0.1796544 
## ... Procrustes: rmse 0.01109283  max resid 0.03652654 
## Run 13 stress 0.218491 
## Run 14 stress 0.4059784 
## Run 15 stress 0.4064148 
## Run 16 stress 0.181569 
## Run 17 stress 0.1847227 
## Run 18 stress 0.2179301 
## Run 19 stress 0.2156235 
## Run 20 stress 0.1815697 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_uuni_during.tif", res = 150, width = 18, 
     height = 7, units = "in")
uunifd
dev.off()
## quartz_off_screen 
##                 2

PERMANOVAs - After samples

Filter samples

Filter after samples and recalculate Bray-Curtis and Jaccard metrics

relab_fecal_filtered_after = relab_fecal_filtered %>% 
  filter(Study_period == "After")
relab_oral_after = relab_oral %>% 
  filter(Study_period == "After")

brayfa = vegdist(relab_fecal_filtered_after[,2:496], method = "bray")
jaccfa = vegdist(relab_fecal_filtered_after[,2:496], method = "jaccard")
brayoa = vegdist(relab_oral_after[,2:496], method = "bray")
jaccoa = vegdist(relab_oral_after[,2:496], method = "jaccard")

Bray-Curtis

Fecal

adonis2(brayfa ~ ethnicity, data = relab_fecal_filtered_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayfa ~ ethnicity, data = relab_fecal_filtered_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1   0.4762 0.07579 2.3782 0.0226 *
## Residual  29   5.8065 0.92421                
## Total     30   6.2827 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(brayoa ~ ethnicity, data = relab_oral_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayoa ~ ethnicity, data = relab_oral_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1  0.24789 0.14292 2.5014 0.0174 *
## Residual  15  1.48655 0.85708                
## Total     16  1.73444 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Bray-Curtis dissimilarity with after samples.

## Run 0 stress 0.1590347 
## Run 1 stress 0.1944615 
## Run 2 stress 0.1937004 
## Run 3 stress 0.1895439 
## Run 4 stress 0.205105 
## Run 5 stress 0.1658979 
## Run 6 stress 0.1713295 
## Run 7 stress 0.1711856 
## Run 8 stress 0.1590348 
## ... Procrustes: rmse 6.513386e-05  max resid 0.0001955462 
## ... Similar to previous best
## Run 9 stress 0.2038524 
## Run 10 stress 0.1714717 
## Run 11 stress 0.1590347 
## ... Procrustes: rmse 0.0001124624  max resid 0.0003000413 
## ... Similar to previous best
## Run 12 stress 0.1658979 
## Run 13 stress 0.1590347 
## ... Procrustes: rmse 9.894018e-05  max resid 0.000277995 
## ... Similar to previous best
## Run 14 stress 0.1866295 
## Run 15 stress 0.2021307 
## Run 16 stress 0.2031112 
## Run 17 stress 0.1656963 
## Run 18 stress 0.1658979 
## Run 19 stress 0.1632431 
## Run 20 stress 0.196653 
## *** Solution reached

## Run 0 stress 0.1367903 
## Run 1 stress 0.1635865 
## Run 2 stress 0.1637332 
## Run 3 stress 0.1367904 
## ... Procrustes: rmse 7.20447e-05  max resid 0.0001608161 
## ... Similar to previous best
## Run 4 stress 0.1504065 
## Run 5 stress 0.1504066 
## Run 6 stress 0.1367903 
## ... Procrustes: rmse 3.209441e-05  max resid 8.068756e-05 
## ... Similar to previous best
## Run 7 stress 0.1617068 
## Run 8 stress 0.1345672 
## ... New best solution
## ... Procrustes: rmse 0.0376222  max resid 0.1126213 
## Run 9 stress 0.1769398 
## Run 10 stress 0.1345672 
## ... Procrustes: rmse 2.563029e-05  max resid 5.312082e-05 
## ... Similar to previous best
## Run 11 stress 0.1367903 
## Run 12 stress 0.1345672 
## ... Procrustes: rmse 1.463215e-05  max resid 3.969615e-05 
## ... Similar to previous best
## Run 13 stress 0.1345672 
## ... Procrustes: rmse 1.981793e-05  max resid 4.80818e-05 
## ... Similar to previous best
## Run 14 stress 0.1519729 
## Run 15 stress 0.1617068 
## Run 16 stress 0.1635865 
## Run 17 stress 0.1617068 
## Run 18 stress 0.1367905 
## Run 19 stress 0.1563232 
## Run 20 stress 0.1520582 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_bray_after.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(brayfa + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(brayfa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(brayoa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Jaccard

Fecal

adonis2(jaccfa ~ ethnicity, data = relab_fecal_filtered_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccfa ~ ethnicity, data = relab_fecal_filtered_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1   0.6328 0.07096 2.2149 0.0156 *
## Residual  29   8.2848 0.92904                
## Total     30   8.9175 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(jaccoa ~ ethnicity, data = relab_oral_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccoa ~ ethnicity, data = relab_oral_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1  0.36814 0.11723 1.9919 0.0196 *
## Residual  15  2.77222 0.88277                
## Total     16  3.14036 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Jaccard similarity with after samples.

## Run 0 stress 0.1590347 
## Run 1 stress 0.1658979 
## Run 2 stress 0.2130876 
## Run 3 stress 0.2118508 
## Run 4 stress 0.2053138 
## Run 5 stress 0.1871108 
## Run 6 stress 0.1712887 
## Run 7 stress 0.1590348 
## ... Procrustes: rmse 0.0001206207  max resid 0.0003348898 
## ... Similar to previous best
## Run 8 stress 0.1714715 
## Run 9 stress 0.1632431 
## Run 10 stress 0.1714068 
## Run 11 stress 0.1926743 
## Run 12 stress 0.1773822 
## Run 13 stress 0.1866586 
## Run 14 stress 0.1590347 
## ... Procrustes: rmse 3.398221e-05  max resid 9.946736e-05 
## ... Similar to previous best
## Run 15 stress 0.1816165 
## Run 16 stress 0.1867049 
## Run 17 stress 0.1590348 
## ... Procrustes: rmse 5.010818e-05  max resid 0.0001729842 
## ... Similar to previous best
## Run 18 stress 0.1656963 
## Run 19 stress 0.1590348 
## ... Procrustes: rmse 0.0001394517  max resid 0.0004449221 
## ... Similar to previous best
## Run 20 stress 0.1912909 
## *** Solution reached

## Run 0 stress 0.1367902 
## Run 1 stress 0.1504065 
## Run 2 stress 0.1617068 
## Run 3 stress 0.1345672 
## ... New best solution
## ... Procrustes: rmse 0.03748586  max resid 0.1125334 
## Run 4 stress 0.1367902 
## Run 5 stress 0.1563232 
## Run 6 stress 0.1563232 
## Run 7 stress 0.1367909 
## Run 8 stress 0.1367905 
## Run 9 stress 0.1637331 
## Run 10 stress 0.1367902 
## Run 11 stress 0.1577926 
## Run 12 stress 0.1637332 
## Run 13 stress 0.1367907 
## Run 14 stress 0.1367904 
## Run 15 stress 0.1769394 
## Run 16 stress 0.1367902 
## Run 17 stress 0.1367904 
## Run 18 stress 0.1563232 
## Run 19 stress 0.1345672 
## ... Procrustes: rmse 2.443653e-05  max resid 4.620881e-05 
## ... Similar to previous best
## Run 20 stress 0.1345672 
## ... New best solution
## ... Procrustes: rmse 2.487242e-05  max resid 4.789656e-05 
## ... Similar to previous best
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_jacc_after.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(jaccfa + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(jaccfa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(jaccoa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Weighted UniFrac

Fecal

adonis2(wuni_fecala ~ ethnicity, data = relab_fecal_filtered_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_fecala ~ ethnicity, data = relab_fecal_filtered_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1  0.09394 0.07012 2.1868 0.0774 .
## Residual  29  1.24577 0.92988                
## Total     30  1.33971 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(wuni_orala ~ ethnicity, data = relab_oral_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_orala ~ ethnicity, data = relab_oral_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.07966 0.11412 1.9323 0.1236
## Residual  15  0.61840 0.88588              
## Total     16  0.69806 1.00000

NMDS visualization of weighted Unifrac distances with after samples.

## Run 0 stress 0.103874 
## Run 1 stress 0.1030343 
## ... New best solution
## ... Procrustes: rmse 0.03049715  max resid 0.1471977 
## Run 2 stress 0.127867 
## Run 3 stress 0.1283418 
## Run 4 stress 0.1340227 
## Run 5 stress 0.1108314 
## Run 6 stress 0.1079067 
## Run 7 stress 0.1079067 
## Run 8 stress 0.1079067 
## Run 9 stress 0.103874 
## Run 10 stress 0.1038767 
## Run 11 stress 0.1254346 
## Run 12 stress 0.1079067 
## Run 13 stress 0.1079068 
## Run 14 stress 0.1079068 
## Run 15 stress 0.1079068 
## Run 16 stress 0.112468 
## Run 17 stress 0.1079067 
## Run 18 stress 0.1328807 
## Run 19 stress 0.103874 
## Run 20 stress 0.1079068 
## Run 21 stress 0.1079067 
## Run 22 stress 0.1038762 
## Run 23 stress 0.1030343 
## ... New best solution
## ... Procrustes: rmse 8.352202e-05  max resid 0.0003853117 
## ... Similar to previous best
## *** Solution reached

## Run 0 stress 0.05926746 
## Run 1 stress 0.07191275 
## Run 2 stress 0.08726035 
## Run 3 stress 0.07191275 
## Run 4 stress 0.05926746 
## ... Procrustes: rmse 1.319554e-06  max resid 4.077547e-06 
## ... Similar to previous best
## Run 5 stress 0.08358197 
## Run 6 stress 0.05926746 
## ... Procrustes: rmse 1.720742e-06  max resid 4.009741e-06 
## ... Similar to previous best
## Run 7 stress 0.05926746 
## ... New best solution
## ... Procrustes: rmse 4.973794e-06  max resid 1.59931e-05 
## ... Similar to previous best
## Run 8 stress 0.07401827 
## Run 9 stress 0.07401828 
## Run 10 stress 0.09924667 
## Run 11 stress 0.07401827 
## Run 12 stress 0.07191275 
## Run 13 stress 0.07191275 
## Run 14 stress 0.07191275 
## Run 15 stress 0.08605518 
## Run 16 stress 0.09687146 
## Run 17 stress 0.09478662 
## Run 18 stress 0.07401827 
## Run 19 stress 0.05926746 
## ... New best solution
## ... Procrustes: rmse 9.301701e-07  max resid 2.343873e-06 
## ... Similar to previous best
## Run 20 stress 0.07191275 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_wuni_after.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(wunifa + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(wunifa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(wunioa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Unweighted UniFrac

Fecal

adonis2(uuni_fecala ~ ethnicity, data = relab_fecal_filtered_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_fecala ~ ethnicity, data = relab_fecal_filtered_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)   
## ethnicity  1  0.20279 0.07361 2.3043 0.0044 **
## Residual  29  2.55215 0.92639                 
## Total     30  2.75494 1.00000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(uuni_orala ~ ethnicity, data = relab_oral_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_orala ~ ethnicity, data = relab_oral_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.04386 0.06782 1.0914 0.3346
## Residual  15  0.60285 0.93218              
## Total     16  0.64671 1.00000

NMDS visualization of unweighted UniFrac distances with after samples.

## Run 0 stress 0.1657305 
## Run 1 stress 0.2285569 
## Run 2 stress 0.1658235 
## ... Procrustes: rmse 0.005360905  max resid 0.02528318 
## Run 3 stress 0.1789747 
## Run 4 stress 0.165719 
## ... New best solution
## ... Procrustes: rmse 0.001621189  max resid 0.005848108 
## ... Similar to previous best
## Run 5 stress 0.2010731 
## Run 6 stress 0.1678215 
## Run 7 stress 0.1824019 
## Run 8 stress 0.182164 
## Run 9 stress 0.1681404 
## Run 10 stress 0.2237202 
## Run 11 stress 0.1675407 
## Run 12 stress 0.1678212 
## Run 13 stress 0.223468 
## Run 14 stress 0.1681404 
## Run 15 stress 0.1781267 
## Run 16 stress 0.1658233 
## ... Procrustes: rmse 0.005827622  max resid 0.02552845 
## Run 17 stress 0.2010731 
## Run 18 stress 0.1657301 
## ... Procrustes: rmse 0.001529493  max resid 0.005840392 
## ... Similar to previous best
## Run 19 stress 0.1904371 
## Run 20 stress 0.1657305 
## ... Procrustes: rmse 0.001633068  max resid 0.005843246 
## ... Similar to previous best
## *** Solution reached

## Run 0 stress 0.1364074 
## Run 1 stress 0.1543538 
## Run 2 stress 0.1367323 
## ... Procrustes: rmse 0.04383167  max resid 0.1382546 
## Run 3 stress 0.1304384 
## ... New best solution
## ... Procrustes: rmse 0.04509471  max resid 0.1268594 
## Run 4 stress 0.1384731 
## Run 5 stress 0.1364076 
## Run 6 stress 0.1385508 
## Run 7 stress 0.1438585 
## Run 8 stress 0.1543538 
## Run 9 stress 0.1304384 
## ... Procrustes: rmse 2.125227e-05  max resid 6.013072e-05 
## ... Similar to previous best
## Run 10 stress 0.1304388 
## ... Procrustes: rmse 0.0009354498  max resid 0.002662176 
## ... Similar to previous best
## Run 11 stress 0.1364078 
## Run 12 stress 0.1364078 
## Run 13 stress 0.3689216 
## Run 14 stress 0.1363188 
## Run 15 stress 0.1304383 
## ... New best solution
## ... Procrustes: rmse 0.000102872  max resid 0.0002883375 
## ... Similar to previous best
## Run 16 stress 0.1694527 
## Run 17 stress 0.1384728 
## Run 18 stress 0.1304383 
## ... New best solution
## ... Procrustes: rmse 0.0003611925  max resid 0.001027821 
## ... Similar to previous best
## Run 19 stress 0.1380133 
## Run 20 stress 0.1367322 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_uuni_after.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(uunifa + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(uunifa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(uunioa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Alpha Diversity

Alpha diversity was calculated using species-level metaphlan results in vegan (2.5-7). Shannon diversity and Shannon-Weaver evenness (Pielou’s evenness) were calculated.

shannond_fecal = diversity(relab_fecal_filtered[,2:496], index = "shannon")
shannond_oral = diversity(relab_oral[,2:496], index = "shannon")

shannone_fecal = shannond_fecal/log(specnumber(relab_fecal_filtered[,2:496]))
shannone_oral = shannond_oral/log(specnumber(relab_oral[,2:496]))

fecal_diversity = cbind(relab_fecal_filtered[,1], relab_fecal_filtered[,497:520], shannond_fecal, shannone_fecal)
fecal_diversity$Study_period = as.factor(fecal_diversity$Study_period)
oral_diversity = cbind(relab_oral[,1], relab_oral[,497:519], shannond_oral, shannone_oral)
oral_diversity$Study_period = as.factor(oral_diversity$Study_period)

A linear model was used to examine differences between self-reported ethnicities and study period (before/during/after) in diversity.

div_fecal = lm(shannond_fecal ~ ethnicity*Study_period, 
                         data = fecal_diversity)
summary(div_fecal)
## 
## Call:
## lm(formula = shannond_fecal ~ ethnicity * Study_period, data = fecal_diversity)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1.3741 -0.2227  0.1292  0.3666  0.6455 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            2.825245   0.117163  24.114   <2e-16 ***
## ethnicityCaucasian                     0.023833   0.180926   0.132    0.895    
## Study_periodBefore                    -0.001497   0.190236  -0.008    0.994    
## Study_periodDuring                     0.117951   0.143495   0.822    0.413    
## ethnicityCaucasian:Study_periodBefore -0.202191   0.293399  -0.689    0.492    
## ethnicityCaucasian:Study_periodDuring  0.002998   0.222444   0.013    0.989    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4971 on 105 degrees of freedom
## Multiple R-squared:  0.03296,    Adjusted R-squared:  -0.01309 
## F-statistic: 0.7157 on 5 and 105 DF,  p-value: 0.613
Anova(div_fecal)
## Anova Table (Type II tests)
## 
## Response: shannond_fecal
##                         Sum Sq  Df F value Pr(>F)
## ethnicity               0.0023   1  0.0094 0.9229
## Study_period            0.7211   2  1.4591 0.2371
## ethnicity:Study_period  0.1599   2  0.3236 0.7243
## Residuals              25.9444 105
summary(glht(div_fecal,linfct=mcp(Study_period="Tukey")))
## Warning in mcp2matrix(model, linfct = linfct): covariate interactions found --
## default contrast might be inappropriate
## 
##   Simultaneous Tests for General Linear Hypotheses
## 
## Multiple Comparisons of Means: Tukey Contrasts
## 
## 
## Fit: lm(formula = shannond_fecal ~ ethnicity * Study_period, data = fecal_diversity)
## 
## Linear Hypotheses:
##                       Estimate Std. Error t value Pr(>|t|)
## Before - After == 0  -0.001497   0.190236  -0.008    1.000
## During - After == 0   0.117951   0.143495   0.822    0.687
## During - Before == 0  0.119449   0.171249   0.698    0.763
## (Adjusted p values reported -- single-step method)
even_fecal = lm(shannone_fecal ~ ethnicity*Study_period, 
                          data = fecal_diversity)
summary(even_fecal)
## 
## Call:
## lm(formula = shannone_fecal ~ ethnicity * Study_period, data = fecal_diversity)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.31337 -0.03586  0.02678  0.07622  0.18643 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            0.649303   0.024716  26.271   <2e-16 ***
## ethnicityCaucasian                     0.005016   0.038167   0.131    0.896    
## Study_periodBefore                     0.030051   0.040131   0.749    0.456    
## Study_periodDuring                     0.019407   0.030270   0.641    0.523    
## ethnicityCaucasian:Study_periodBefore -0.086399   0.061893  -1.396    0.166    
## ethnicityCaucasian:Study_periodDuring -0.003113   0.046925  -0.066    0.947    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1049 on 105 degrees of freedom
## Multiple R-squared:  0.03584,    Adjusted R-squared:  -0.01007 
## F-statistic: 0.7807 on 5 and 105 DF,  p-value: 0.5658
Anova(even_fecal)
## Anova Table (Type II tests)
## 
## Response: shannone_fecal
##                         Sum Sq  Df F value Pr(>F)
## ethnicity              0.00359   1  0.3266 0.5689
## Study_period           0.01186   2  0.5395 0.5846
## ethnicity:Study_period 0.02733   2  1.2427 0.2928
## Residuals              1.15454 105
summary(glht(even_fecal,linfct=mcp(Study_period="Tukey")))
## Warning in mcp2matrix(model, linfct = linfct): covariate interactions found --
## default contrast might be inappropriate
## 
##   Simultaneous Tests for General Linear Hypotheses
## 
## Multiple Comparisons of Means: Tukey Contrasts
## 
## 
## Fit: lm(formula = shannone_fecal ~ ethnicity * Study_period, data = fecal_diversity)
## 
## Linear Hypotheses:
##                      Estimate Std. Error t value Pr(>|t|)
## Before - After == 0   0.03005    0.04013   0.749    0.732
## During - After == 0   0.01941    0.03027   0.641    0.796
## During - Before == 0 -0.01064    0.03613  -0.295    0.953
## (Adjusted p values reported -- single-step method)
div_oral = lm(shannond_oral ~ ethnicity*Study_period, data = oral_diversity)
summary(div_oral)
## 
## Call:
## lm(formula = shannond_oral ~ ethnicity * Study_period, data = oral_diversity)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.32262 -0.10843 -0.04301  0.11353  0.66514 
## 
## Coefficients:
##                                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            3.37870    0.06839  49.403   <2e-16 ***
## ethnicityCaucasian                    -0.05358    0.10658  -0.503    0.619    
## Study_periodBefore                    -0.04956    0.09672  -0.512    0.612    
## ethnicityCaucasian:Study_periodBefore  0.06253    0.15073   0.415    0.681    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2163 on 30 degrees of freedom
## Multiple R-squared:  0.01195,    Adjusted R-squared:  -0.08685 
## F-statistic: 0.121 on 3 and 30 DF,  p-value: 0.9471
Anova(div_oral)
## Anova Table (Type II tests)
## 
## Response: shannond_oral
##                         Sum Sq Df F value Pr(>F)
## ethnicity              0.00410  1  0.0877 0.7691
## Study_period           0.00482  1  0.1030 0.7504
## ethnicity:Study_period 0.00805  1  0.1721 0.6812
## Residuals              1.40317 30
even_oral = lm(shannone_oral ~ ethnicity*Study_period, data = oral_diversity)
summary(even_oral)
## 
## Call:
## lm(formula = shannone_oral ~ ethnicity * Study_period, data = oral_diversity)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.062819 -0.016022 -0.001432  0.020062  0.057567 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            0.724772   0.009848  73.599   <2e-16 ***
## ethnicityCaucasian                    -0.030771   0.015346  -2.005    0.054 .  
## Study_periodBefore                    -0.008127   0.013926  -0.584    0.564    
## ethnicityCaucasian:Study_periodBefore  0.031961   0.021703   1.473    0.151    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03114 on 30 degrees of freedom
## Multiple R-squared:  0.124,  Adjusted R-squared:  0.03645 
## F-statistic: 1.416 on 3 and 30 DF,  p-value: 0.2574
Anova(even_oral)
## Anova Table (Type II tests)
## 
## Response: shannone_oral
##                           Sum Sq Df F value Pr(>F)
## ethnicity              0.0018014  1  1.8576 0.1830
## Study_period           0.0002153  1  0.2221 0.6409
## ethnicity:Study_period 0.0021031  1  2.1687 0.1513
## Residuals              0.0290921 30

Shannon diversity and evenness did not vary across self-reported ethnicities or study periods in either fecal or oral samples.

## quartz_off_screen 
##                 2
## quartz_off_screen 
##                 2